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A  general  scattering  formulation  is  presented  for  predicting  the  far-field  scat¬ 
tered  pressure  from  irregular,  axisymmetric,  finite-length  bodies  for  three  bound¬ 
ary  conditions — soft,  rigid  and  fluid.  The  formulation  is  an  extension  of  a  two- 
dimensional  conformal  mapping  approach  [D.  T.  DiPerna  and  T.  K.  Stanton,  J. 
Acoust.  Soc.  Am.  96:  3064-3079  (1994)]  to  scattering  by  finite-length  bodies. 
This  extended  formulation,  which  is  inherently  numerically  efficient  to  evaluate,  in¬ 
volves  conformally  mapping  the  surface  of  an  irregular,  finite-length  body  to  a  new, 
orthogonal  coordinate  system  in  which  the  separation  of  variables  method  may  be 
used  to  solve  the  Helmholtz  equation  and  satisfy  the  boundary  conditions.  Extensive 
comparisons  with  previously  published  results  using  other  formulations  are  presented. 
This  formulation  is  shown  to  be  very  accurate  in  the  prediction  of  scattering  from 
smooth,  symmetric  bodies  for  a  wide  range  of  frequencies  (Rayleigh  through  geomet¬ 
ric  scattering  region),  scattering  angles  (monostatic  and  bistatic),  aspect  ratios,  and 
for  each  of  the  three  boundary  conditions  listed  above.  Reasonable  agreement  has 
also  been  demonstrated  for  irregular,  realistic  shapes  with  soft  boundary  conditions. 
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I.  INTRODUCTION 


The  prediction  of  acoustic  scattering  from  finite  and  infinitely  long  bodies  has  been  pursued 
for  many  years,  starting  with  Lord  Rayleigh’s  work  on  scattering  from  a  sphere  (Rayleigh,  1945). 
Exact  analytical  solutions  to  the  acoustic  wave  equation  require  the  scatterer’s  surface  to  exactly 
match  the  locus  of  all  points  for  which  the  radial  coordinate  is  a  constant.  Such  exact  analytical 
solutions  exist  only  for  a  limited  number  of  cases  for  which  the  separation  of  variables  is  possible 
(Morse  and  Feshbach,  1953;  Bowman  et  a/.,  1987).  In  all  of  these  cases,  the  boundary  is  simple; 
e.g.,  a  sphere,  infinitely  long  cylinder  and  prolate  spheroid. 

For  complex  shapes,  approximate  analytical  solutions,  including  the  perturbation  method, 
and  approximate  asymptotic  formulations,  such  as  physical  optics  (Gaunaurd,  1985)  and  the 
geometric  theory  of  diffraction  (Levy  and  Keller,  1959;  Yamashita,  1990)  have  been  developed. 
Numerical  solutions  have  also  been  developed,  including  the  boundary  element  method  (Tobac- 
man,  1984;  FYancis,  1993),  T-matrix  (Waterman,  1968;  Varadan  et  a/.,  1982;  Lakhtakia  et  a/., 
1984;  Hackman  and  Todoroff,  1985)  and  the  mode  matching  methods  (Yamashita,  1990).  All  of 
these  approaches  are  limited  in  one  or  more  of  the  following:  frequency  range,  class  of  surfaces, 
types  of  boundary  conditions,  eccentricity  of  shape  and/or  computational  implementation  and 
numerical  efficiency. 

DiPerna  and  Stanton  (1994)  introduced  a  conformal  mapping  approach  to  predicting  far- 
field  sound  scattering  by  infinitely  long  cylinders  of  noncircular  cross  section.  The  approach, 
termed  the  Fourier  Matching  Method  (FMM),  involves  a  conformal  mapping  of  variables  to  a 
new  coordinate  system  in  which  the  constant  radial  coordinate  exactly  matches  the  scatterer 
surface.  The  method  makes  use  of  the  Newton-Raphson  algorithm  to  execute  the  mapping. 
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The  boundary  conditions  are  satisfied  by  requiring  the  Fourier  coefficients  in  the  new  angular 
variable  of  the  total  field  to  be  zero,  and  then  the  resultant  scattered  field  is  expressed  in  terms 
of  circular  eigenfunctions. 

The  FMM  proved  to  be  accurate  over  a  wide  range  of  frequencies,  shapes  of  cross  section, 
and  penetrable  (fluid)  as  well  as  impenetrable  boundary  conditions.  Furthermore,  the  approach 
is  inherently  numerically  efficient  due  to  the  nature  of  its  formulation.  For  example,  the  FMM 
was  shown  by  DiPerna  and  Stanton  (1994)  to  be  more  efficient  than  the  T-matrix  method  for  the 
case  of  the  high-aspect-ratio  elliptic  cylinder  because  fewer  terms  were  needed  for  the  numerical 
integrations.  Even  after  incorporating  the  FMM  basis  functions  into  the  T-matrix  calculations, 
the  FMM  required  85%  fewer  integration  points. 

A  major  limitation  of  the  two-dimensional  FMM  was  the  fact  that  it  was  formulated  for 
the  case  of  an  infinitely  long  cylinder— a  two-dimensional  scattering  solution.  Many  practical 
scattering  problems  involve  scattering  from  finite  bodies  and  cannot  be  accurately  modeled  by 
the  two-dimensional  solution.  In  order  to  address  this  need,  the  FMM  is  extended  in  this  pa¬ 
per  to  predict  the  scattering  from  finite-length  bodies.  In  order  for  this  particular  approach 
to  be  used  for  finite  bodies,  the  outer  boundary  of  the  bodies  must  be  described  by  a  function 
rotated  about  the  length-wise  axis.  Hence,  although  the  function  is  arbitrary  and  these  bodies 
are  three-dimensional,  they  are  restricted  to  axisymmetric  shapes.  As  with  the  two-dimensional 
formulation,  this  approach  is  intrinsically  numerically  efficient  and  is  valid  over  a  wide  range  of 
frequencies  and  shapes  as  well  as  both  monostatic  and  bistatic  scattering  geometries.  The  ex¬ 
tension  has  been  formulated  for  three  boundary  conditions — Dirichlet  (soft,  or  pressure- release) , 
Neumann  (rigid)  and  Cauchy  (fluid).  In  all  three  cases,  the  surrounding  material  is  fluid. 
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In  Section  II,  the  theoretical  basis  for  the  formulation  is  presented,  which  includes  the  devel¬ 
opment  of  the  new  orthogonal  coordinate  system  to  which  the  body  is  mapped,  the  conformal 
mapping  procedure,  modal  series  solutions  to  the  transformed  Helmholtz  equation,  and  resulting 
equations  for  the  modal  series  coefficients  after  satisfying  the  three  boundary  conditions.  In 
Section  III,  a  description  of  the  numerical  implementation  of  the  FMM  is  presented,  as  well  as 
several  practical  numerical  issues  that  arise  in  the  solution  of  the  scattering  problem.  In  Section 
IV,  the  numerical  results  are  presented  for  various  shapes  (spheres,  smooth  prolate  spheroids  and 
two  finite  bodies  with  irregular  surfaces),  boundary  conditions  (soft,  rigid  and  fluid),  and  over  a 
wide  range  of  frequencies  and  scattering  angles.  The  results  axe  compared  with  various  previ¬ 
ously  published  results  using  other  approaches.  Section  V  contains  a  summary  and  concluding 
remarks. 

II.  THEORY 

The  derivation  of  the  extended  formulation  for  scattering  by  an  axisymmetric  finite-length 
body  is  conceptually  very  similar  to  the  corresponding  derivation  of  the  two-dimensional  solution 
described  by  DiPerna  and  Stanton  (1994);  in  fact,  some  of  the  elements  are  identical.  Both 
solutions  begin  with  the  wave  equation  in  a  known  coordinate  system  and  conformally  map  the 
coordinate  variables  to  a  new,  orthogonal  coordinate  system  in  which  the  locus  of  all  points 
where  the  new  radial  coordinate  is  a  constant  exactly  coincides  with  the  scatterer  surface.  The 
difference  in  the  coordinate  systems  between  the  two  cases  concerns  the  fact  that  one  involves 
two-dimensional  coordinates  while  the  other  involves  an  additional  coordinate  dimension  with 
a  new  geometry  defined  for  the  finite  body.  Both  solutions  use  mapping  functions  that  are 
identical  in  form,  transform  the  Helmholtz  equation  to  the  new  coordinate  system,  and  satisfy 
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the  boundary  conditions  using  identical  techniques  to  arrive  at  differing,  yet  structurally  similar, 
expressions  for  the  scattered  pressure.  The  two-dimensional  solution  includes  a  mapping  function 
that  corresponds  to  the  shape  of  the  boundary  of  a  cross-sectional  slice  of  the  cylinder,  while 
the  three-dimensional  solution  uses  a  mapping  function  that  corresponds  to  the  shape  of  the 
boundary  of  a  length-wise  slice  of  the  body  (specifically,  the  function  that  is  rotated  about  the 
longitudinal  axis).  Furthermore,  the  two-dimensional  solution  includes  circular  eigenfunctions 
while  the  three-dimensional  solution  for  the  scattered  pressure  is  expressed  in  terms  of  spherical 
wave  functions;  i.e.,  spherical  Bessel  and  Hankel  functions  and  associated  Legendre  functions. 
Due  to  these  similarities,  the  original  work  will  be  referred  to  quite  regularly  in  the  development 
that  follows. 

Consider  the  scalar  wave  equation: 


V2P  = 


1  d2P 
c2  dt2  ’ 


(1) 


where  P  is  the  acoustic  pressure  in  three  dimensions,  V2  is  the  Laplacian  operator,  c  is  the  speed 
of  sound,  and  t  is  time.  Assuming  a  harmonic  time  dependence,  e-1"*,  where  u  is  the  angular 
frequency,  the  wave  equation  becomes  the  scalar  Helmholtz  differential  equation  in  Cartesian 
coordinates: 

V2P(®,  y ,  z)  +  k2P(x,  y,  z)  -  0.  (2) 

Here,  k  =  w/c  =  2-k/X  is  the  spatially  independent  acoustic  wave  number,  and  A  is  the  acoustic 
wavelength.  In  all  cases  considered  here,  the  body  does  not  support  shear  waves  and  is  sur¬ 
rounded  by  a  fluid  medium.  By  a  conformal  (angle-  and  orientation-preserving)  transformation 
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of  coordinates,  the  transformed  Helmholtz  equation  in  the  new  coordinate  system  becomes: 


V2P(u,  w ,  v)  +  k2F(u,  w)P(uy  w,  v)  =  0,  (3) 

where  (u,  w ,  v)  are  the  new  coordinates,  and  F(u}  w)  is  a  function  which  depends  on  the  specific 
transformation  (Morse  and  Feshbach,  1953;  DiPerna  and  Stanton,  1994).  With  the  exception 
that  the  wave  number  is  now  a  function  of  position,  the  new  Helmholtz  equation  is  formally 
identical  to  the  Helmholtz  equation  in  Cartesian  coordinates. 

A.  Conformal  mapping 

Since  x,  y  and  z  are  mutually  orthogonal  in  the  Cartesian  coordinate  system,  conformally 
mapping  them  into  a  new  coordinate  system  guarantees  that  the  new  coordinates  (ii,tu,v)  will 
be  mutually  orthogonal,  which  eases  the  computation  of  the  normal  particle  velocity  on  the 
boundary.  Additionally,  the  conformal  mapping  generates  a  new  set  of  angular  functions  which 
fit  the  scatterer  surface  more  naturally;  that  is,  points  along  the  surface  that  change  rapidly 
in  (x,y,*)  are  sampled  at  a  higher  spatial  rate  yet  are  equally  spaced  in  (u,  w,  t/).  A  new 
coordinate  system  must  first  be  established,  and  then  the  conformal  mapping  function  is  defined 
and  expanded  to  provide  a  method  by  which  the  body  may  be  mapped  to  the  new  coordinate 
system. 

1.  Orthogonal  Coordinate  System 

An  orthogonal  coordinate  system  can  be  generated  for  a  three-dimensional  body  of  revolution 
from  a  two-dimensional  conformal  mapping.  Consider  the  geometry  in  Fig.  1,  in  which  <f>  is  the 
azimuthal  angular  coordinate  ranging  from  0  to  2 it  (measured  from  the  positive  x-axis  in  the 
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xy-plane),  6  is  the  polax  angular  coordinate  ranging  from  0  to  v  (measured  from  the  positive 
z-axis),  and  r  is  the  radial  coordinate  ranging  from'O  to  oo.  This  body  is  one  of  revolution  that 
is  formed  by  rotating  the  contour  of  the  body  about  the  z  axis,  in  the  same  way  that  the  prolate 
spheroidal  coordinate  system  is  created  from  an  ellipse  rotated  about  the  major  axis  (Flammer, 
1957).  Consider  a  new  coordinate  system  whose  azimuthal  angular  coordinate,  v,  corresponds 
to  <f>  in  the  original  coordinate  system.  The  new  polar  angular  coordinate,  to,  is  measured  from 
the  polar  axis,  z,  and  ranges  from  0  to  tt,  as  does  the  original  polar  angular  coordinate,  6.  The 
scatterer  surface  in  the  original  coordinate  system  is  defined  by  the  vector,  ~r  ,  but  in  the  new 
coordinate  system  the  scatterer  surface  is  defined  by  the  locus  of  all  points  where  the  new  radial 
coordinate  is  a  constant;  specifically,  u  =  0.  Defining  the  functions,  f(u,  to)  and  g(u,  to),  of  the 
new  coordinate  system  as  shown  in  Fig.  1,  trigonometry  prescribes  dimensions  of  the  body  in 
the  x,  y  and  z  directions  to  be: 


x(u,  to,  v)  =  f(u,  to)  cos(v) 

(4) 

y(u,  to,  v )  =  f(u,  to)  sin(v) 

(5) 

z(u,w,v)  =  g(u,w). 

(6) 

The  position  vector,  ~r,  is  defined  in  the  new  coordinate  system  by: 

7*  ( u ,  to,  v)  =  x(u,  to,  v)i  4-  y(u,  to,  v)j  +  z(ti,  to,  v)k,  (7) 
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where  t,  j,  and  k  are  unit  vectors  along  the  coordinate  axes.  The  position  vector  can  be 
alternatively  expressed  by  substituting  Eqs.  (4)-(6)  into  Eq.  (7): 

T*  =  /(u,  ui)  cos(v)i  +  f(u,  w)  sin  (v)j  +  g(u,  w)k.  (8) 

The  local  projection  of  Ir  in  each  of  the  coordinate  directions  is  given  by  the  partial  derivative 
of  ~r  with  respect  to  each  of  the  variables: 

=  /u(u,  to)  cos(v)?  +  /U(u,  to)  Bin(v)J + gu(u,  w)k  (9) 

~?xu  =  fw(u,  w)<x>s(v)7+ fw(u,w)sm(v)j+gw(u,w)k  (10) 

=  -/(«,  to)  sin(«)t  +  f(u,  to)  cos(v)J,  (11) 

where  the  subscript  denotes  the  variable  with  respect  to  which  the  partial  derivative  is  taken. 

As  mentioned  earlier ,  an  orthogonal  coordinate  system  is  desirable  since  it  facilitates  the  com¬ 
putation  of  the  normal  particle  velocity  on  the  boundary  necessary  for  satisfying  the  boundary 
conditions.  An  orthogonal  coordinate  system  requires  the  following  condition  to  be  satisfied: 

^•^  =  0  (12) 

(13) 

* 

(14) 


T  tU  *  T  y  *—  0 

r\,  •  r\„  =  0, 
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which  can  be  expanded  as: 


/U(u,  w)f(u,  w)  cos(v)  sin(u)(-l  +  1)  =  0  (15) 

fw(u,w)f(u,  w)  cos(v)  sin(u)(— 1  +  1)  =  0  (16) 

fu(u,  w)  fw{u, w)(cos2(v)  +  sm2(v))  +  gu(u,w)gw(u,w)  =  0.  (17) 

The  first  two  conditions  are  automatically  satisfied.  The  third  condition  simplifies  to: 

/u(u,tu)/w(u,te)  +  gu(u,w)gv](u,w)  =  0,  (18) 

which  will  be  satisfied  if: 

fu(u,w)  =  gw(u,w)  (19) 

and 

U{u,v))  = -gu(u,  to).  •  (20) 

These  are  precisely  the  Cauchy- Riemann  equations  for  an  analytic  function  (Hildebrand,  1964). 
Therefore,  if  f(u, «)  and  g(u,  v )  are  chosen  to  be  harmonic,  then  the  Cauchy- Riemann  conditions 
will  be  satisfied,  making  them  analytic  functions  which  represent  a  conformed  transformation. 
A  shape  initially  plotted  in  the  (x,  y,  z)  coordinate  system  will  be  transformed  into  a  shape  in 
the  (u,  w,v)  coordinate  system  with  changes  in  position  and  size  while  preserving  angles  and 
proportions  (Morse  and  Feshbach,  1953).  Orthogonality  of  the  coordinate  system  as  well  as  the 
form  of  the  Helmholtz  equation  will  be  preserved  (Strang,  1986). 


2.  Mapping  Function 

As  discussed  above,  a  conformal  mapping  function  must  be  developed  to  map  the  scatterer 
from  the  old  coordinate  system  to  the  new  orthogonal,  axisymmetric  coordinate  system  just 
developed.  It  must  be  noted  at  this  point  that,  to  the  authors’  knowledge,  a  general  three- 
dimensional  mapping  does  not  exist  in  the  field  of  mathematics.  Due  to  the  fact  that  conformal 
mappings  are  currently  limited  to  two  dimensions,  the  geometry  for  the  finite  body  must  be 
axisymmetric  about  one  of  the  axes.  The  particular  mapping  used  herein  is  a  two-dimensional 
mapping  developed  by  DiPema  and  Stanton  (1994)  extended  to  a  finite  body  of  revolution  which 
is  axisymmetric  about  the  longitudinal  axis.  The  infinitely  long  cylindrical  geometry  in  DiPema 
and  Stanton  (1994)  was  described  in  circular  cylindrical  coordinates  with  the  radial  coordinate, 
r,  being  a  function  of  6,  the  azimuthal  angular  coordinate  ranging  from  0  to  27r.  The  conformal 
mapping  in  that  case  applied  to  the  function,  r(6),  which  corresponded  to  the  shape  of  the 
boundary  of  a  cross-sectional  slice.  In  this  work,  6  is  now  the  polar  angular  coordinate  ranging 
from  0  to  7r,  and  4>  is  the  azimuthal  angular  coordinate  ranging  from  0  to  2i r.  The  function, 
r(0),  and  associated  conformal  mapping  is  now  associated  with  the  shape  of  the  boundary  of  a 
length-xoise  slice. 

The  category  of  surfaces  described  by  Eqs.  (4)-(6)  has  the  additional  limitation  that  r  be 
single- valued;  i.e.,  there  can  be  only  one  value  of  r  for  each  w.  Following  DiPema  and  Stanton 
(1994),  the  mapping  procedure  for  the  axisymmetric  finite  body  is  commenced  by  expanding  r 
in  a  Fourier  series  relative  to  the  polar  angle,  9,  shown  in  Fig.  1: 

00 

r(6)  =  a  +  ]T[r£  cos  (rtf)  -I-  r*n  sin(n0)],  (21) 

n=l 
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where  a  is  the  average  radius  of  the  body,  and  r*  and  r*  are  the  usual  Fourier  series  coefficients 
that  in  this  case  correspond  to  the  deviation  of  the  surface  from  the  shape  of  a  circle.  Note 
that  the  series  requires  more  terms  to  converge  for  a  high  aspect  ratio  (ratio  of  length  to  width) 
prolate  spheroid  compared  to  the  Fourier  series  for  a  shape  that  varies  little  from  the  shape  of 
a  circle.  Rewriting  the  cos(n0)  and  sin(»0)  functions  in  terms  of  exponentials  and  using  the 
expression: 

=  ^  [**«  +  *<).  (22) 

gives 

|r|  =  aei6  +  f^lK^1+n)e  +  (23) 

«=1 

For  a  conformal  mapping  from  the  (®,  y,  z)  coordinate  system  to  the  new  coordinate  system 
in  (tt,u>,t>): 

M(p)  =  M(u  +  iw),  (24) 

where  M(p)  is  the  analytic  mapping  function  in  terms  of  it,  the  radial  variable,  and  to,  the 
polar  angular  variable,  and  p  =  u  +  iw.  It  is  desirable  to  make  scattering  predictions  using 
this  model  without  inversely  mapping  the  results  of  this  model  back  to  the  original  coordinate 
system.  The  potentially  difficult  inverse  mapping  is  avoided  by  choosing  M(p)  such  that  the 
coordinate  system  becomes  spherical  as  the  radial  coordinate  is  increased.  While  the  choice 
of  such  a  mapping  function  allows  predictions  of  this  model  to  be  easily  compared  to  existing 
solutions,  it  restricts  direct  comparisons  to  the  far-field  only.  The  general  approach  can  certainly 
be  used  in  the  near-field,  but  comparisons  of  near-field  scattering  between  this  formulation  and 
other  solutions  would  require  an  inverse  mapping.  Note  also  that  there  will  be  two  different 
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mapping  functions:  G(p)  for  the  exterior  problem  and  T(p)  for  the  interior  problem.  Fbr  the 
interior  problem,  T(p)  is  chosen  such  that  the  coordinate  system  becomes  spherical  as  the  radial 
coordinate  is  decreased. 

For  the  exterior  problem,  G{p)  must  be  chosen  such  that: 

(1)  As  u  — +  oo,  the  coordinate  system  becomes  spherical, 

(2)  the  transformed  Helmholtz  equation  is  solvable,  and 

(3)  u  =  0  is  the  scatterer  surface. 

The  first  two  conditions  can  be  satisfied  by  choosing  the  form  of  the  exterior  mapping  function 
(DiPema  and  Stanton,  1994)  to  be: 


00 

G(p)  =  c_,e'  +  ]Tcne-’*,  (25) 

n=0 

which  can  be  decomposed  into  the  complex  components: 

oo 

<?(u,  w)  s  Re(G(p))  =  c_ieu  cos(u»)  4-  c„e~nu  cos(mu)  (26) 

n=0 

and 

OO 

/(u,  w)  =  Im(G(p))  =  c_ie“  Bin(tz>)  +  ^  Cne~nu  sin(nti;).  (27) 

n=0 

The  coefficients  of  the  mapping  function  must  be  chosen  such  that  u  =  0  defines  the  scatterer 
surface  in  the  new  coordinate  system.  Equating  the  functions  in  Eq.  (23)  and  Eq.  (25)  (with 
u  =  0)  that  represent  the  surface  in  the  two  coordinate  systems  gives: 

oo  oo 

oe"  +  ^[/£ei(1+n>9  +  =  c_1e’w  +  £  cne-’'nu’.  (28) 

»=1  n*»0 
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Since  the  left-hand  side  contains  positive  and  negative  frequency  components  while  the  right- 
hand  side  contains  only  negative  frequency  components  (with  the  exception  of  c_  ie,w),  0  and  w 
are  not  equal;  therefore,  it  is  necessary  to  determine  the  extent  to  which  6  depends  on  w.  Since 
it  was  assumed  earlier  that  the  surface  is  periodic  and  can  be  represented  as  a  Fourier  series,  the 
deviation  of  9  from  w  will  be  periodic  and  can  be  represented  as  a  Fourier  series.  Specifically, 
assume: 

00 

6{w)  =  tu  +  cos(Ztu)  +  6f  sin(Zio)] .  (29) 

'=i 

The  conformal  mapping  relies  on  the  choice  of  Sf  and  6*  such  that  Eq.  (28)  is  satisfied.  Us¬ 
ing  the  orthogonality  relationships  of  complex  exponential  functions,  multiplying  both  sides  by 
(l/2?r)e_<,t"  and  integrating  over  w  from  0  to  2tt  gives: 

i.  jf e-Hv>  dw  =  j  <  j  } '  (30) 

where  j  is  an  integer.  This  set  of  nonlinear  constraints  is  identical  in  form  to  that  of  DiPerna 
and  Stanton  (1994)  and  is  solved  by  use  of  an  extension  of  the  Newton-Raphson  method,  the 
details  of  which  are  laid  out  in  Appendix  A  of  DiPerna  and  Stanton  (1994).  Note,  however,  that 
even  though  the  integral  in  Eq.  (30)  is  performed  from  0  to  2?r,  tv  is  defined  in  the  scattering 
geometry  from  0  to  ir  only  (not  2ir).  Consequently,  the  mapping  coefficients  are  computed  based 
on  the  periodic  extension  from  0  to  2ir,  but  only  half  of  them  are  used.  The  upper  result  in  the 
right  hand  side  of  Eq.  (30)  is  used  to  solve  for  the  values  of  Sf  and  Sf,  which  are  then  used  to 
solve  for  the  mapping  coefficients,  Cn,  through  use  of  the  lower  result  in  the  right  hand  side  of 
Eq.  (30). 
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The  uniqueness  of  the  transformation  is  tested  by  verifying  that  the  Jacobian  of  the  transfor¬ 
mation  is  nonzero.  This  ensures  that  there  exists  only  one  (x,z)  for  each  (u,w).  Specifically: 

|G'(p)|2  ¥=  0,  u  >  0.  (31) 

The  interior  mapping  procedure  is  identical  to  the  exterior  mapping  procedure  with  the 
exception  that  the  interior  mapping  function,  T(p)y  is  different  from  G{p)\  specifically,  T(p)  is 
chosen  such  that  the  coordinate  system  becomes  spherical  as  the  radial  coordinate  is  decreased. 
The  remaining  two  conditions  mentioned  above  in  the  choice  of  G(p)  remain  the  same  for  the 
interior  problem.  These  conditions  are  satisfied  by: 

oo 

T(p)  =  Y,tnenp.  (32) 

n—O 

To  summarize,  the  procedure  described  above  to  conformally  map  the  scatterer  shape  from 
the  original  coordinate  system  to  a  new  coordinate  system  is  identical  in  form  to  the  procedure 
presented  in  DiPerna  and  Stanton  (1994).  In  this  study,  the  same  mapping  procedure  is  extended 
to  a  different  (finite-length,  axisymmetric)  scattering  geometry;  specifically,  it  is  extended  to  the 
shape  of  the  boundary  in  the  length-wise  slice.  The  results  of  this  mapping  will  be  used  in 
solving  the  Helmholtz  equation  in  three  dimensions  in  the  next  section. 

B.  Solutions  to  the  Helmholtz  equation 

The  3-dimensional  Helmholtz  equation  from  Eq.  (2)  in  the  original  coordinate  system  is 
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expressed  in  spherical  coordinates  as: 


V2P(r,M)  +  *2^M^)  =  0, 


(33) 


the  general  solution  to  which  is: 

P“‘(r>M)  =  X]  ]C  anmJn(^)Pnm(cos(0))e<m^  (34) 

n=— oo  m=— oo 

+  E  E 

n=— oo  m=— oo 

where  jn(kr)  is  the  spherical  Bessel  function  of  the  first  kind  of  order  n,  h£\kr)  is  the  spherical 
Hankel  function  of  the  first  kind  of  order  n,  and  P^(cos(9))  is  the  associated  Legendre  function 
of  degree  n  and  order  m.  The  radial  coordinate  is  r,  the  polar  angular  coordinate  is  0,  and 
the  azimuthal  angular  coordinate  is  <f>.  The  scattered  field  coefficients,  6nm,  are  determined  by 
satisfying  the  boundary  conditions  using  the  known  coefficients,  Onm,  of  the  incident  plane  wave 
field  traveling  from  the  6q  direction  relative  to  the  2-axis: 

Onm  -  *”€m(2n  +  1)  ||p^(cos(g0)),  (35) 

where  em  is  the  Neumann  factor,  and  V  is  the  gamma  function. 

P“*(r,0,0)  is  the  total  pressure  external  to  the  scatterer  in  the  original  coordinate  system 
(i.e.  before  transformation):  the  first  term  in  Eq.  (34)  represents  the  incident  pressure,  and  the 
second  term  represents  the  scattered  pressure.  Quantities  in  the  original  coordinate  system  can 
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be  expressed  in  terms  of  the  new  coordinate  system  defined  in  Sec.  II.A.l  and  Fig.  1: 


Pat(u,w,v)  =  £  £  anmMkriu,™))^  e<mU 

n^oom^oo  \r(u,w)  J 

+  E  E  e‘mv, 

n=-oom=-oo  if;;  y 


which  is  now  the  expression  for  the  total  far-held  pressure  in  the  new  coordinate  system  as  the 
sum  of  the  incident  and  scattered  pressure  fields,  respectively. 

The  procedure  to  determine  the  pressure  field  inside  the  scatterer  is  identical  in  nature  to 
the  exterior  problem,  but  the  mapping  function  is  different  for  the  interior  problem;  therefore, 
the  new  Helmholtz  equation  is  identical  in  form,  but  incorporates  a  different  wave  number,  kit 
to  accurately  characterize  the  material  properties  of  the  scatterer’s  interior.  The  expression  for 
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the  interior  pressure  field  becomes: 


00  oo 

Pin\u,W,v)  =  ^  5Z  lnmjn{k\r{u,  Vj))Pn 

n=— oo  m=— oo 


(41) 


where  /nm  axe  the  internal  field  coefficients.  Only  the  spherical  Bessel  function  is  included  in 
the  expression  for  the  internal  field  since  the  spherical  Hankel  function  becomes  infinite  at  the 
origin. 

In  the  limit  of  great  distances  from  the  scatterer,  the  coordinate  system  becomes  spherical, 
and  the  asymptotic  form  of  the  Hankel  function  varies  inversely  with  distance.  The  scattered 
pressure  (the  second  term  in  Eq.  (40))  in  this  fax-field  limit  has  the  form: 

pikr 

pscat  (42) 

u— *oo  T 


where  the  scattering  amplitude,  fs,  having  units  of  length,  is  a  measure  of  the  efficiency  with 
which  an  object  scatters  sound  and  is  a  function  of  the  object’s  size,  shape,  orientation,  material 
properties  and  the  wavelength  of  the  incident  wave.  The  scattering  amplitude  is  given  in  general 
form  as: 


oo  00 

/.-EE  w— 

n=— oo  m=-oo 


eimv. 


(43) 


The  far-field  scattered  energy  evaluated  in  the  backscatter  direction  is  often  expressed  in 
terms  of  the  target  strength  (TS)  with  units  of  decibels  (dB)  relative  to  1  m  (Urick,  1983),  which 
is  given  by: 

TS  =  10  log  au  =  10  log  \fbtf ,  (44) 
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where  /**  is  the  scattering  amplitude  evaluated  in  the  backscattering  direction,  and  ab,  =  |/t*|2 
is  the  differential  backscattering  cross  section,  which  differs  from  the  often-used  backscattering 
cross  section,  a,  by  a  factor  of  An  ( a  =  Ancrbt).  In  order  to  compare  scattering  from  objects  of 
different  sizes  but  similar  proportions,  target  strength  is  often  normalized  by  the  the  square  of 
some  typical  dimension.  Using  the  length  ( L )  of  the  elongated  scatterer  as  the  normalization 
constant,  the  “reduced”  target  strength  (RTS)  is: 


RTS  =  10  log  (2j)  =  lOlogl/4,12  -  10  log(L2)  =  10  log 


(45) 


In  the  case  of  a  sphere,  the  target  strength  is  normalized  by  7 ra2  instead  of  L2.  An  alternative 
expression  often  used  to  represent  the  energy  scattered  in  three  dimensions  is  the  normalized, 
steady-state  form  function  (Neubauer,  1986): 


(46) 

where  a  is  the  radius  of  a  sphere  or  the  semi-minor  axis  of  a  prolate  spheroid. 

To  summarize,  a  new  orthogonal  coordinate  system  has  been  established  to  which  the  original 
coordinate  system  is  mapped  via  a  conformal  mapping  function.  The  scatterer  surface  in  the 
new  coordinate  system  is  defined  by  the  locus  of  all  points  where  the  radial  coordinate  is  a 
constant.  The  solution  to  the  Helmholtz  equation  in  the  new  coordinate  system  will  have  the 
same  form  as  the  original  solution,  but  with  new  coordinates  using  the  realizations  of  Eqs.  (36) 
(38).  This  can  be  confirmed  by  inspecting  the  case  of  the  exterior  field  for  the  spherical 
scatterer:  r(u,to)  =  c^eu>  0(u,w)  =  w  and  <j>  =  t;,  resulting  in  the  same  solution  as  in  spherical 
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coordinates  using  separation  of  variables.  Also,  for  the  case  far  (u  — >  oo)  from  the  non-spherical 
scatterer:  r(u,w)  -*  c_ie“,  0(u,w)  -*  w  and  <f>  =  v,  resulting  in  the  same  solution  as  in 
spherical  coordinates  in  the  far-field.  Equations  (40)  and  (41)  are  a  general  solution  for  the 
total  pressure  in  the  case  of  acoustic  scattering  from  a  finite-length,  axisymmetric  body  for  all 
frequencies,  all  angles  (bistatic  scattering),  and  for  impenetrable  (soft  and  rigid)  and  penetrable 
(fluid)  boundary  conditions.  All  the  quantities  needed  to  satisfy  soft,  hard  and  fluid  boundary 
conditions  have  now  been  established.  Each  boundary  condition  generates  a  different  set  of 
modal  series  coefficients,  bnm,  which  are  used  in  Eq.  (43)  to  compute  the  scattering  amplitude. 
These  coefficients  will  be  determined  for  three  boundary  conditions  in  the  following  section. 

C.  Boundary  conditions 

The  solution  for  the  scattered  field  in  Eqs.  (40)  and  (43)  depends  on  the  boundary  conditions. 
For  each  boundary  condition,  there  is  a  different  set  of  scattered  field  coefficients,  bnm,  which 
needs  to  be  evaluated.  These  coefficients  are  determined  in  the  normal  approach  in  which  the 
pressures  and/or  velocities  are  matched  at  the  boundaries.  This  matching  is  specific  to  the 
material  properties  on  each  side  of  the  boundary.  For  each  of  the  different  boundary  conditions, 
expressions  are  derived  in  the  following  sections  for  &nm»  all  of  which  are  in  terms  of  the  parameters 
of  the  new  coordinate  system. 

The  challenge  presented  by  an  irregular  body  is  finding  a  solution  to  a  system  of  equations  in 
order  to  satisfy  the  boundary  conditions;  specifically,  a  set  of  functions,  ip,  must  be  determined 
by  which  the  basis  functions  are  multiplied  to  generate  the  system  of  equations.  For  a  separable 
geometry,  ip  is  chosen  to  be  the  angular  eigenfunction  to  yield  a  closed-form  solution  based  on 
the  orthogonality  of  the  eigenfunctions.  In  the  case  of  irregular  surfaces  that  do  not  conform 
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to  coordinate  surfaces,  a  set  of  functions,  tp)  must  be  chosen  to  solve  the  system  of  equations, 
DiPerna  and  Stanton  (1994)  chose  the  eigenfunctions  in  the  new  (cylindrical)  coordinate  system, 
~  e~tmv.  It  is  in  this  choice  of  eigenfunctions  that  the  FMM  derives  its  name.  Without 
the  conformal  change  of  variables  to  (uiw)v)}  this  choice  of  functions  would  not  be  possible. 
For  the  three-dimensional  case,  the  eigenfunctions  in  the  spherical  coordinate  system  are  chosen; 
specifically  (Morse  and  Feshbach,  1953): 


C  =  P?(  cos  W)  sinMe-*™.  (47) 

Using  this  choice  for  a  system  of  equations  is  generated  to  satisfy  each  of  the  boundary 
conditions  discussed  herein  (soft,  rigid  and  fluid).  The  system  of  equations  is  solved  for  the 
scattered  field  coefficients,  6nm,  which  are  then  used  in  Eq.  (43)  to  determine  the  scattering 
amplitude. 


1.  Dirichlet  boundary  conditions 

For  soft  (pressure  release  or  Dirichlet)  boundary  conditions,  the  total  pressure  vanishes  on  the 
surface  of  the  scatterer;  i.e.,  =  0.  Using  Eq.  (40),  the  series  solution  for  the  total 

exterior  pressure  field  evaluated  at  the  boundary  with  pressure  release  boundary  conditions  is 
set  equal  to  zero: 


n^ocmtr^  \r(U0,w)J 

+  Y,  £  W?)(Mu0,tfl))^f^)eiTO  =  0. 

n=-oo  m=-oo  V  W)  J 


(48) 
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The  system  of  equations  necessary  to  satisfy  this  boundary  condition  is  generated  by  multiplying 
both  sides  of  this  equation  by  from  Eq.  (47)  and  integrating  over  the  range  of  w  and  v  (Morse 
and  Feshbach,  1953): 


+  t  t 

n--ooms“00  x  v  ' '  / 

x  P™  (cos  (w))  sin  (w)  e_,mu}  dvdw  =  0. 


(49) 


Performing  the  integration  on  v  gives: 

°nm  (jn(kr(uQ,  w))P™  P^(cos(u;))sm(^)^  dw  (50) 

+  £  6""*/  i7(C08(w))8lll(w))  dw  =  0, 

for  each  azimuthal  order  m.  This  system  of  equations  for  pressure  release  boundary  conditions 
can  be  written  in  compact  form: 

oo  oo 

£  a„m/C+  £  bnmQ™  —  0,  (51) 

oo  «=- oo 

where  a„m  is  given  in  Eq.  (35),  and  and  Q™  are  defined,  using  Eq.  (50),  as: 

V  =  jT  "))V  Jr(co»W)  *(•))  dw  (52) 
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«r = [  c(cosw)sinw)  <*«.. 


(53) 


Using  Eq.  (51),  the  series  coefficients,  bnm>  for  the  scattered  field  due  to  a  soft  boundary  are: 

bnm  —  ~  {Q™)  *  IC^m]  (54) 

where  ( . , .  .)_1  denotes  a  matrix  inversion. 

2.  Neumann  boundary  conditions 

With  rigid,  or  Neumann,  boundary  conditions,  the  normal  particle  velocity  vanishes  on  the 
scatterer  surface;  i.e.,  n  •  VPext{uo,w,v)  =  0,  where  n  is  the  unit  vector  normal  to  the  surface, 
and  V  is  the  gradient  operator.  The  method  to  solve  for  bnm  is  broadly  similar  to  the  previous 
case,  although  several  new  terms  are  necessary  to  satisfy  this  particular  boundary  condition. 
The  unit  vectors  in  the  u,  «,  and  w  directions  are  (Hildebrand,  1964): 


if  -  r“ 

<lu  —  pqr“T 

!*•«! 

(55) 

— ► 

^  —  '  tu 

(Lui  —  T 

\rw\ 

(56) 

— ► 

&v  —  | — ►  |  • 

I**. 

(57) 

Using  Eqs.  (9)  -  (11),  the  scale  factors  of  each  of  the  new  coordinates  (Morse  and  Feshbach, 
1953)  can  be  rewritten  using  the  new  coordinates: 

hu  -  |r\,|  =  V72(u, w)+gl(u,w)  =  y/ (58) 
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K  =  |rg  =  y/ f£{%u>)  +  glfaw)  =  \/ fl{u,w)  + 


K  =  |rM  =  /( u,w). 


Note  that: 


hu(u,w )  =  hw(u,w)  =  h(u,w). 


In  the  new  coordinate  system,  the  unit  vectors  become: 


8.  =  £d^cos(v)i  +  M±4sm(v)3+ 

h(u,w )  v  1  h(u,w)  "  h(u,  w) 


S.  =  #24cosW?+  #!4Sm(»)?+  #^4* 
fc(u,w)  v  '  "  h(u,w) 


ay  =  —  sin(u)t  4-  cos(v)j. 


The  scatterer  surface  has  been  defined  by  u  =  uq\  therefore,  n  =  o^.  The  gradient  of  P  is: 


nn/  ,  1  dP„  1  dP^  ,  1  dP 

{ u,w,v  h^faw)  du0™  hw(u,w)dwaw  hviUfW)  dv 


Using  n  =  Eu  and  /it,(u,u>)  =  1/ /£(u,  to)  +  f£{u,  w): 


n  •  VP(u,i»,«)  = 


«0  +  />.  *«)  5u 


Using  Eq.  (40),  the  series  solution  for  the  normal  particle  velocity  of  the  total  external  field 


for  rigid  boundary  conditions  evaluated  at  the  boundary  is  set  equal  to  zero: 


s' v  ( (ra)  <-) 

+B' v  {Lij-hSHkrM)r?  (ra)  H  -  °- 


The  system  of  equations  necessary  to  satisfy  this  boundary  condition  is  generated  by  multiplying 
both  sides  of  the  equation  by  ij)™  in  Eq  (47)  and  integrating  over  the  range  of  w  and  v: 

xPn{cos(w))6in(w)e~imv}  dvdw  =  0. 

Using  Eq.  (66)  and  performing  the  integration  on  v  gives: 

OO  y 

E"™./  (jn(kr(uo,w))Pnu  ( ^4)  ^  W)  -  9(M0,  uK(uo,  w) 

r^oo  Jo  \  r»(tto,w)  (69) 

+  £  Km  r  (hW{kr(uo,w))P™  ( w)9u(vo, w)  - g(up, w)ru(uo, w) 

Jo  \  n'u  \r(uo,w)J  r*K«/) 


where,  as  with  the  other  terms,  the  subscript,  u,  also  indicates  the  partial  derivative  of  the 
Legendre  functions  with  respect  to  u. 
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As  in  the  previous  section,  the  system  of  equations  for  rigid  boundary  conditions  can  be 


written  in  compact  form: 


£  <wC+  £  bnmQT  =  0, 


where  <w  is  given  in  Eq.  (35),  R!?  is  the  integral  in  the  first  half  of  Eq.  (69),  and  Q?  is  the 
integral  in  the  second  half  of  Eq.  (69).  Using  Eq.  (70),  the  series  coefficients,  bnmf  for  the 
far-field  scattered  field  due  to  a  rigid  boundary  are: 


3.  Cauchy  boundary  conditions 

Fluid  (Cauchy)  boundary  conditions  is  the  case  in  which  the  body  consists  of  a  fluid  with 
material  properties  other  than  that  of  the  surrounding  fluid  and  does  not  support  a  shear  wave. 
In  this  case,  the  exterior  and  interior  pressure  fields  and  normal  components  of  particle  velocity 
are  required  to  be  equal  on  the  surface.  The  fluid  boundary  conditions  take  the  form: 

Pressure: 


t  t  W.(M^))a”(^)e'~ 

n=— oo  m=- oo  N  v  /  7 

n=—  oo  m=s— oo  ' 

-£  E 

n=-oo  m=-oo  x  7 


25 


Particle  Velocity: 


'  IS,  i>^->)*(£$)*-) 

+B- v  (£$)  e‘m') 

= s- v  C£JL  w^*.»«r  fej) '-) . 


(73) 


where  lnm  are  the  interior  field  coefficients.  The  system  of  equations  necessary  to  satisfy  this 
boundary  condition  is  generated  by  multiplying  both  sides  of  the  equation  by  VC  in  Eq  (47)  and 
integrating  over  the  range  of  w  and  v: 

Pressure: 


($$) e*"") 

xPC(  cosH)  sin(io)e_,mv}  dvdw 

i  r  { ($$) 

xP?( oob(w))  sin(w)e-<mu}  dvdw 


(74) 
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Particle  Velocity: 


LT{h  (£3) ' 

t  ($$)<-)) 

\n=-oom=-*oo  x  x  /7  77 


x  P™(cos(w))  sin(iu)e  ,mu }  dvdu> 


a:  {bv  (i^'""j"(Mu°'w))c  (k^o 


xPnm(cos(u»))  sin(^o)e-,m,'}  dvdw. 


Using  Eq.  (66)  and  performing  the  integration  on  v  for  both  Eqs.  (74)  and  (75)  gives: 

Pressure: 

E  jn(kr(uo,w))P?  P™{cas(w))sm(w)dw 

+  flbnmJ0  hn)(kr(uo^VJ))PT  (r(^C))  (76) 

=  E  j»(Muo.™))7T  P?(coa(w))8in(w)dw 
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Particle  Velocity: 


E  «-  /"Yj.(*r(«o. «))?.  rKl»)i.(„l,)-a(%»),,K.) 

»=-oo  V  ’  \r(uo>w)/  r2(u0,«>) 


:  sin(io)dio 


+A,(M«o.«))^  fcr„(u0ltu)')  fy(co3(ui))  sin(u,),ftu 

+  f'bmf( AW(trK,  -)?  “)  -  ■ “fr-fo-  ■")  (77) 

*=-«,  -/o  \  ’  r(uo,u>)  r2(uo,ty)  v  ’ 


($5$)  *.K»))  ^jjjS 


i^*(cos(u»)) 


<")  +  /2Kw) 


sw(w)dw 


n=-^oo  ^  V  ’  \r(uo,w)J  r2(uo,w) 

+*-(*.>•  («..«-)«*  *•.(«.,»)) 

Vr(tio,w);  /  ^Kto)  +  /2(uo,w)  v  ' 

Writing  Eqs.  (76)  and  (77)  in  compact  form: 

00  00  OO 

53  ttnm^n  +  53  ^  tim'SJT  (78) 

n=-oo  n=— oo  n— — oo 

iC  a™n^»n  53  ^nm^n*  “  53  ^nm<®nn »  (79) 

**-oo  n=— oo  n=— oo 

where  is  given  in  Eq.  (35);  J2*  Q™  and  are  now  the  first,  second  and  third  integrals  in 
Eq.  (76);  and  iJJJ*1 ,  Q%*  and  ££”  are  now  the  first,  second  and  third  integrals  in  Eq.  (77).  In 
contrast  to  the  above  two  cases  of  impenetrable  boundary  conditions,  bnm  must  be  solved  for  in 
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two  simultaneous  equations.  Solving  for  bnm  in  Eqs.  (78)  and  (79): 

bnm = -  (Qr-s?  -  Qnsrr1  (ks?  -  jc$d<w  (so) 

In  summary,  the  scattering  amplitude  (Eq.  (43))  of  a  particular  body  is  determined  by 
conformally  mapping  the  scatterer  surface  to  a  new  coordinate  system  by  solving  Eq.  (30),  then 
solving  for  the  scattering  coefficients,  bnmi  for  the  appropriate  boundary  condition  as  specified 
in  the  equations  of  this  section. 

III.  NUMERICAL  IMPLEMENTATION 

A.  General  approach 

Prediction  of  the  scattering  by  an  axisymmetric  finite-length  body  using  the  above  formulation 
requires  a  series  of  steps  involving  numerical  methods.  Given  the  complexity  of  the  procedure, 
it  is  summarized  briefly  in  this  section,  with  some  of  the  challenges  elaborated  upon  in  the  next 
section.  The  shape  of  the  scatterer  is  first  described  by  the  array  of  points,  (r,0),  in  polar 
coordinates  in  the  zz-plane  (Fig.  1).  The  terms  r  and  6  are  then  expanded  in  a  Fourier  series 
using  Eqs.  (21)  and  (29).  The  conformal  mapping  is  based  upon  determining  the  values  of  6^ 
and  6*,  which  are  determined  by  solving  the  upper  result  in  the  right  hand  side  of  Eq.  (30). 
This  is  a  set  of  nonlinear  constraints  which  are  solved  in  an  iterative,  numerical  manner  by  an 
extension  of  the  Newton- Raphson  method.  The  details  of  the  Newton- Raphson  method  can 
be  found  in  the  Appendix  of  DiPerna  and  Stanton  (1994)  or  in  many  math  texts.  Once  the 
values  of  6f  and  6*  are  determined,  the  lower  result  in  the  right  hand  side  of  Eq.  (30)  is  used  to 
determine  the  mapping  coefficients,  c„,  which  are  then  used  to  compute  functions,  g(u,  w)  and 
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/( u,tw),  of  the  new  coordinate  system  (Eqs.  (26)  and  (27)).  After  the  values  of  g(u,w)  and 

/(u,  m)  are  determined,  the  conformal  mapping  is  complete,  and  the  solution  to  the  Helmholtz 

equation  in  the  new  coordinate  system  can  be  determined.  The  number  of  terms  included  in 

the  summations  in  Eqs.  (21),  (29)  and  (30)  are  chosen  in  an  iterative  manner — predictions  are 

made  at  first  with  a  small  number  of  terms,  and  then  the  number  is  increased  for  subsequent 

predictions  until  the  scatterer  is  accurately  mapped  to  the  new  coordinate  system.  The  mapping 

can  be  verified  visually  by  using  g{uy  w)  and  /(u,  it/)  to  plot  on  top  of  the  shape  in  the  original 

coordinate  system  (Fig.  1).  Fbr  some  shapes,  the  conformal  mapping  method  can  be  avoided 

by  using  a  function  to  describe  the  boundary.  In  the  case  of  a  smooth  prolate  spheroid,  the 

functions  /  =  osin(u/)  and  g  =  6cos(iu)  can  be  used  in  Eqs.  (26)  and  (27),  where  a  and  b  are 

the  semi- minor  and  semi-major  axes  of  the  prolate  spheroid,  respectively.  In  the  new  coordinate 

system,  the  external  far-field  pressure  (Eq.  (40))  is  solved  by  use  of  Eqs.  (36)  -  (38).  For 

soft  boundary  conditions,  solve  Eq.  (54)  by  numerical  integration  of  Eqs.  (£2)  and  (53).  The 

numerical  integration  can  be  performed  by  use  of  one  of  numerous  integration  techniques,  many 

of  which  are  described  in  Press  et  ai.,  (1992).  In  this  work,  simple  matrix  summation  was 

used.  Likewise  for  rigid  and  fluid  boundary  conditions,  solve  Eqs.  (71)  and  (80)  by  performing 

the  numerical  integrations  in  Eqs.  (69)  and  (76)  -  (77),  respectively.  With  the  scattered  field 

,‘n 

coefficients,  bnmy  the  scattering  amplitude  is  finally  computed  Eq.  (43).  The  point  at  which  the 

A 

summations  in  Eq.  (43)  are  truncated  is  determined  in  an  iterative  manner.  The  issues  and 
challenges  involved  in  truncation  and  performing  the  matrix  inversions  in  Eqs.  (54),  (71)  and 
(80)  are  discussed  in  the  next  section. 
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G.  Numerical  issues 


In  general,  there  exist  significant  inherent  difficulties  in  the  numerical  implementation  of  in¬ 
finite  series  solutions.  For  series  solutions  for  simple  geometries  such  as  a  sphere  or  cylinder, 
and  with  currently  available  computers,  the  series  converges  relatively  rapidly.  However,  as  the 
shape  deviates  from  a  simple  geometry,  higher  modes  are  required  for  a  converged  solution.  The 
problem  is  particularly  complicated  for  objects  of  high  aspect  ratio  and  irregularity  (roughness) 
on  the  boundary.  The  accuracy  of  the  solution  must  be  balanced  with  accounting  for  the  degree 
of  roughness  and  elongation  of  the  boundary. 

1.  Precision 

Some  of  the  difficulty  in  the  numerical  calculation  of  the  wave  functions  is  merely  the  time 
required  for  the  computer  to  compute  the  wave  functions.  Beyond  this  issue  of  time,  the  limiting 
factors  in  the  numerical  implementation  are  machine  precision  and  matrix  manipulation. 

The  FMM  generates  a  transition  matrix,  much  like  the  T-matrix  model  (Waterman,  1968), 
that  relates  the  incident  field  coefficients  to  the  scattered  field  coefficients.  For  a  spherical 
scatterer,  the  transition  matrix  is  diagonal  and  each  non-zero  term  on  the  main  diagonal  is 
an  eigenvalue  for  each  mode  computed.  If  the  scatterer  shape  deviates  from  spherical,  the 
matrix  contains  off-diagonal  terms.  The  additional  higher  modal  terms  required  to  represent 
the  scattering  become  extremely  small,  sometimes  falling  below  the  value  that  can  be  accurately 
represented  numerically,  resul^pg  in  a  singular  matrix  in  which  the  true  values  of  its  elements 
are  below  the  precision  of  the:  machine.  Thus,  machine  roundoff  error  is  introduced  into  the 
solution  and  quickly  dominates  the  results  as  it  propagates  through  the  solution  via  repetitive 
matrix  manipulation. 


31 


The  fundamental  problem  is  finite  machine  precision.  The  number  of  modes  required  to 
accurately  represent  the  scattering  involves  matrix  elements  that  are  sometimes  numerically 
smaller  than  the  machine  can  accurately  compute.  Little  can  be  done  to  improve  precision  since 
it  is  a  hardware  limitation.  Future  technology  may  afford  greater  precision,  but  the  extent  to 
which  it  improves  the  solution  is  difficult  to  estimate.  The  relationship  between  greater  precision 
and  increased  accuracy  in  the  solution  is  not  necessarily  linear,  but  could  be  a  rapidly  decaying 
one  instead. 

Besides  hardware  limitations,  numerical  accuracy  is  very  dependent  upon  accurate  and  effi¬ 
cient  numerical  implementation  of  the  theory.  There  are  a  number  of  numerical  algorithms  that 
improve  matrix  manipulation  and  handling  of  roundoff  error  that  have  been  investigated  in  this 
work,  including  orthogonal  triangular  decomposition,  LU  factorization,  balancing  of  matrices  and 
scaling  strategies  (Press  et  ai}  1992).  The  single  most  important  algorithm  investigated  and 
used  in  the  implementation  of  this  formulation  is  the  well-known  singular  value  decomposition 
(SVD)  algorithm  (Press  et  a/.,  1992). 

For  an  ill-conditioned  matrix  in  which  some  elements  are  below  machine  precision,  the  SVD 
algorithm  sets  those  numerically  indiscernible  elements  to  zero.  Eliminating  one  or  more  linear 
combinations  of  the  set  of  equations  that  is  to  be  solved  is  justified  because  those  subspaces 
are  dominated  by  roundoff  error  and  contribute  negligible  energy  to  the  solution  and  contribute 

principally  to  error;  indeed,  the  error  can  be  amplified.  Singular  values  whose  ratio  to  the 

* 

largest  singular  value  is  less  than  N  times  the  machine  precision  are  set  to  zero  (Press  et  a/., 
1992).  This  threshold  value  is  equivalent  to  the  rank  of  the  matrix,  which  is  an  estimate  of  the 
number  of  linearly  independent  rows  or  columns  of  a  matrix.  Eliminating  erroneous  subspaces 
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yields  a  more  stable  numerical  result  and  a  convergent  solution  at  higher  frequencies;  however, 
a  disadvantage  of  the  method  is  the  potential  of  removing  some  energy  that  contributes  to  the 
solution,  thereby  diminishing  the  amplitude  and  structure  of  the  results  at  higher  frequencies. 
Therefore,  the  threshold  used  for  SVD  must  be  set  with  care. 

2.  Truncation 

Regardless  of  whether  the  wave  functions  are  difficult  to  compute  and  are  susceptible  to  lim¬ 
itations  of  machine  precision,  numerical  computation  of  an  infinite  series  such  as  in  Eqs.  (40) 
and  (43),  cannot  be  performed  exactly  due  to  the  fact  that  there  are  an  infinite  number  of  terms. 
Furthermore,  the  point  at  which  the  series  must  be  truncated  is  not  necessarily  obvious.  For 
lower  frequencies  and  smooth  spherical  or  smooth  low-aspect-ratio  spheroids  (i.e.,  low  eccentric¬ 
ity),  reaching  a  converged  solution  before  the  onset  of  singularity  is  typically  not  a  problem  using 
currently  available  personal  computers.  As  the  frequency,  aspect  ratios  or  degree  of  irregularity 
are  increased,  roundoff  error  begins  to  increase  as  a  factor  in  preventing  the  solution  from  being 
converged. 

Given  the  importance  of  these  numerical  issues,  a  study  was  conducted  to  explore  the  condi¬ 
tions  under  which  the  solution  would  be  reliable  once  numerically  evaluated.  The  “performance 
envelope”  of  the  FMM  for  broadside  backscattering  by  smooth  prolate  spheroids  of  varying  aspect 
ratios  for  soft  and  rigid  boundary  conditions  was  investigated  (Figs.  2  and  3). 

In  this  study,  a  “converged”  solution  is  defined  as  one  in  which  the  computation  of  addi¬ 
tional  modes  does  not  significantly  change  the  result  for  a  given  value  of  ka.  Specifically,  the 
scattering  amplitude  for  a  converged  solution  (although  strictly  a  truncated  form  of  the  exact 
solution)  changed  by  less  than  0.1%  (an  amount  chosen  arbitrarily  in  this  study)  with  the  com- 
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putation  of  additional  modes,  which  correlated  very  well  with  visual  inspection.  “Truncated” 
and  “numerically  stable”  approximations  are  less  objective  and  are  a  relative  indication  of  the 
degree  of  reliability.  Truncated  approximations,  as  defined  in  this  paper,  employ  a  sufficient 
number  of  modes  to  represent  the  scattering  to  a  lesser  degree  of  accuracy  than  the  above-defined 
converged  solution,  but  still  generally  predicts  the  overall  scattering  levels  and  finer  structure. 
Numerically  stable  approximations  use  a  sufficient  number  of  modes  to  adequately  represent  the 
overall  amplitude,  if  not  the  finer  structure,  of  the  majority  of  the  scattering.  For  both  the  soft 
(Fig.  2)  and  rigid  (Fig.  3)  cases,  the  solution  was  obtained  for  prolate  spheroids  with  aspect 
ratios  ranging  from  1:1  (i.e.,  a  sphere)  to  13:1,  and  ka  ranging  from  0  to  10.  Note  that  a  is  the 
semi-minor  axis  of  the  prolate  spheroid,  b  is  the  semi-major  axis,  and  the  aspect  ratio  is  fc/a. 
This  assignment  of  notation  (a  and  6)  differs  from  that  in  some  literature,  as  a  in  this  paper  is 
intended  to  correspond  to  a  cylindrical  radius  for  the  elongated  bodies.  For  each  of  the  two 
boundary  conditions,  a  converged  solution  is  reached  for  the  sphere  and  2:1  aspect  ratio  prolate 
spheroid  for  values  of  ka  up  to  10.  As  the  aspect  ratio  is  increased  beyond  2:1,  the  value  of 
ka  at  which  the  solution  is  converged  falls  off  rapidly  down  to  levels  of  1.3  (soft)  and  1.2  (rigid) 
at  an  aspect  ratio  of  13:1.  The  truncated  approximation  and  numerically  stable  approximation 
behave  less  predictably  as  a  function  of  ka  and  aspect  ratio,  but  indicate  generally  decreasing 
convergence  values  of  ka  with  increasing  aspect  ratio.  The  value  of  ka  at  which  a  converged 
solution  is  reached  in  the  soft  case  is  generally  higher  than  that  of  the  rigid  case,  particularly  for 
aspect  ratios  of  8: 1  or  less.  This  is  consistent  with  the  fact  that  the  soft  case  does  not  include 
an  expression  for  the  derivative  of  pressure  as  in  the  rigid  case  (Eq.  (66)),  thus  requiring  fewer 
matrix  manipulations  that  would  propagate  roundoff  error. 
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IV.  NUMERICAL  RESULTS 

All  of  the  results  plotted  in  Sec.  IV  were  generated  using  Eqs.  (43)-(46),  (54),  (71)  and  (80) 
to  illustrate  predictions  as  well  as  demonstrate  the  accuracy  of  the  FMM  in  a  wide  variety  of 
applications.  Predictions  using  the  FMM  are  first  compared  to  those  using  the  exact  solution 
for  soft,  rigid  and  fluid  spheres.  Predictions  are  then  compared  between  the  FMM  and  various 
formulations  for  smooth  prolate  spheroids  for  a  variety  of  aspect  ratios,  frequencies,  and  incidence 
and  scattering  angles.  Thirdly,  computations  are  presented  for  scattering  from  a  realistic, 
irregular  body.  Lastly,  predictions  are  presented  for  scattering  by  gaseous  smooth  prolate 
spheroids  in  the  region  of  resonance  for  various  aspect  ratios. 

A.  Spheres:  comparison  with  exact  solution 

Computing  the  scattering  from  the  sphere  is  an  essential  benchmark  for  the  FMM  since  the 
exact  solutions  exist.  The  results  for  the  soft,  rigid  and  fluid  spheres  produced  by  the  FMM 
are  identical  to  the  exact  solutions  (Anderson,  1950)  in  all  three  cases  (Fig.  4).  In  the  case  of 
fluid  boundary  conditions,  material  properties  for  a  weak  scatterer  (i.e.,  one  having  properties 
resembling  that  of  the  surrounding  fluid)  were  used.  The  specific  values  were  chosen  to  resemble 
zooplankton  tissue. 

B.  Prolate  spheroids:  comparison  with  various  solutions 

As  mentioned  earlier,  accurately  computing  the  scattering  from  a  prolate  spheroid  is  a  difficult 
task,  particularly  at  high  aspect  ratios  and  high  frequencies.  Every  scattering  model  has  its 
strengths  and  weaknesses,  so  to  thoroughly  test  the  scope  of  applicability  for  the  FMM,  it  is 
compared  to  a  variety  of  models  within  their  ranges  of  performance.  The  approaches  used 
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for  comparison  are  the  Deformed  Finite  Cylinder  (DFC)  model,  the  T-matrix  formulation,  the 
Boundary  Element  Method  (BEM),  the  exact  prolate  spheroidal  solution  and  the  Kirchhoff 
approximation. 

1.  Deformed  finite  cylinder  model  (DFC) 

Stanton  adapted  his  formulation  for  the  finite  cylinder  (Stanton,  1988a,  b)  to  finite  cylinders 
of  deformation;  e.g.,  prolate  spheroids  and  uniformly  bent  finite  cylinders  (Stanton,  1989  ).  The 
model,  which  is  based  on  the  modal  series  solution  for  an  infinitely  long  cylinder,  has  proven 
to  be  accurate  at  broadside  and  near-broadside  orientations  to  the  prolate  spheroid  for  aspect 
ratios  of  5:1  and  higher  (Partridge  and  Smith,  1995).  The  FMM  compares  well  to  the  DFC  for 
the  case  of  broadside  backscatter  from  a  rigid  prolate  spheroid  of  aspect  ratios  2:1,  5:1  and  10:1 
(Fig.  5).  As  the  aspect  ratio  is  increased,  results  for  the  FMM  are  limited  to  lower  frequencies 
to  maintain  a  converged  solution.  In  the  Rayleigh  region  ( ka  <  1),  the  TS  increases  as  expected, 
proportional  to  (ka)4.  Note  also  that  the  agreement  between  the  FMM  and  the  DFC  in  the 
Rayleigh  region  improves  as  the  aspect  ratio  increases,  reflecting  the  fact  that  the  accuracy  of 
the  DFC  improves  with  aspect  ratio,  eventually  coinciding  with  the  FMM  results. 

2.  T-matrix  method 

The  T-matrix  method  is  a  formally  exact,  numerical  solution  that  has  been  the  focus  of  much 
research  to  date  and  has  proven  to  be  an  accurate  model  for  numerous  scattering  problems.  In 
two  earlier  studies,  the  T-matrix  approach  has  been  applied  to  backscatter  and  bistatic  scattering 
by  rigid  prolate  spheroids  with  an  aspect  ratio  of  2:1  (Varadan  et  al.,  1982)  and  10:1  (Hackman, 
1993)  and  will  be  compared  with  FMM  predictions.  There  is  excellent  agreement  between  the 
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FMM  and  the  T-matrix  method  in  the  cases  of  backscatter  for  a  2:1  aspect  ratio  prolate  spheroid 
at  broadside  ( 60  =  90°)  and  oblique  (0O  =  45°)  angles  of  incidence,  as  well  as  in  the  bistatic 
geometry  of  broadside  incidence  ($o  =  90°)  and  end-on  reception  (6r  =  0°)  (Fig.  6).  There  is 
also  excellent  agreement  between  the  T-matrix  and  the  FMM  after  increasing  the  aspect  ratio 
from  2:1  to  10:1  for  the  rigid  prolate  spheroid  in  the  case  of  broadside  and  end-on  backscatter 
(Fig.  7).  All  of  the  cases  in  Figs.  6  and  7  illustrate  the  interference  between  the  specular 
reflection  and  a  Franz,  or  creeping  wave,  and  the  associated  shifts  in  the  peaks  and  nulls  for  the 
different  angles  of  incidence  and  reception  (Uberall  et  al. ,  1966). 

3.  Boundary  element  method  (BEM) 

Francis  (1993)  developed  a  numerical  method  of  computing  the  scattering  by  a  finite  body 
using  a  boundary  element  method  based  on  a  partial  application  of  a  Helmholtz  gradient  formu¬ 
lation.  The  FMM  and  BEM  agree  perfectly  to  high  frequencies  (ka  =  10)  in  the  case  of  end-on 
backscatter  for  a  2:1  aspect  ratio  rigid  prolate  spheroid  (Fig.  8). 

4.  Exact  prolate  spheroidal  solution 

The  exact  solution  for  a  prolate  spheroid  in  prolate  spheroidal  coordinates  has  been  used  by 
various  investigators  (Yeh,  1967;  Furusawa,  1988;  and  Ye  et  al.,  1997).  While  this  is  formally  an 
exact  solution,  it  encounters  the  same  numerical  challenges  mentioned  in  Sec.  III.  Comparisons 
of  the  FMM  to  the  exact  solution  in  the  case  of  backscatter  and  forward  scatter  for  a  range  of 
incidence  angles  (6q  =  30°,  70°,  90°)  demonstrate  good  to  excellent  agreement  in  all  cases  (Fig. 
9). 
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5.  Kirchhoff  approximation 

The  Kirchhoff  approximation  is  a  frequently  used  approximation  for  backscattering;  however, 
it  is  generally  limited  to  the  geometric  scattering  region  (fca  >>  1),  perfectly  reflecting  boundary 
conditions  (soft  and  rigid)  and  near-broadside  incidence  (Born  and  Wolf,  1999).  It  is  a  conve¬ 
nient  approximation  for  the  amplitude  of  the  backscatter,  but  does  not  accurately  represent  the 
oscillations  in  the  amplitude  due  to  the  interaction  between  the  specular  reflections  and  Franz 
waves.  Comparisons  of  the  FMM  to  this  approximation  within  the  range  of  reliability  of  the 
Kirchhoff  approximation  further  establish  the  consistency  of  the  FMM.  Specifically,  the  reduced 
target  strengths  for  broadside  backscattering  generally  agree  well  for  the  single  frequencies  pre¬ 
sented,  with  improved  agreement  as  the  frequency  increases  (Fig.  10).  At  ka  =  1,  the  Kirchhoff 
approximation  is  near  the  limit  of  its  range  of  validity  and  under-estimates  the  amplitude.  Ad¬ 
ditionally,  it  can  be  seen  that  the  Kirchhoff  approximation  falls  off  much  more  rapidly  than  the 
FMM  for  all  fca  as  the  incidence  angle  moves  away  from  broadside. 

C.  Irregular  bodies:  comparison  with  Kirchhoff  approximation 

As  emphasized  earlier,  accurate  computation  of  the  scattering  by  smooth  elongated  bodies 
is  not  a  simple  task,  particularly  as  the  aspect  ratio  and  frequencies  increase.  The  additional 
complexity  of  an  irregular  surface  further  complicates  am  already  difficult  problem.  A  key 
distinction  of  the  FMM  is  its  ability  to  conformally  map  an  irregular  axisymmetric  surface  to  an 
orthogonal  coordinate  system  that  better  fits  the  scatterer  surface.  To  demonstrate  the  practical 
application  of  the  FMM  to  a  realistic,  asymmetrical,  irregular  body  with  non-Gaussian  roughness, 
the  acoustic  scattering  by  a  gas-filled  swimb ladder  from  an  alewife  fish  ( Alosa  pseudoharengus ) 
has  been  modeled  and  compared  with  the  Kirchhoff  approximation  for  two  representative  shapes 
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corresponding  to  the  dorsal  and  ventral  aspects  of  the  swimbladder  (Figs.  11  and  12).  This 
particular  problem  was  chosen  as  fish  are  a  significant  scatterer  of  sound  in  the  ocean,  and  the 
swimbladder  is  the  organ  in  the  fish  that  generally  dominates  the  scattering.  A  more  extensive 
treatment  of  this  application  is  given  in  Reeder  et  al.  (submitted).  At  these  values  of  ka,  soft 
boundary  conditions  were  used,  given  the  high  contrast  in  material  properties  between  the  gas 
and  surrounding  tissue.  There  is  generally  reasonable  agreement  between  FMM  and  Kirchhoff  for 
most  values  of  ka  for  Shape  #1  in  Fig.  11,  but  for  Shape  #2  in  Fig.  12,  the  Kirchhoff  predictions 
generally  under-estimate  the  scattering  relative  to  the  FMM  at  angles  well  off  broadside. 

D.  Low  ka  resonance  scattering  for  gaseous  bodies:  comparison  with  T-matrix  and 
exact  prolate  spheroidal  solutions 

In  contrast  to  the  cases  of  rigid  and  weakly-scattering  objects  described  above  whose  scattering 
levels  decrease  dramatically  in  the  ka  «  1  region,  gaseous  bodies  have  a  strong  resonance  in  that 
region  and  have  been  the  subject  of  many  earlier  studies.  Although  in  the  low  ka  region  there 
are  fewer  restrictions  in  the  calculations,  challenges  still  remain  in  the  prediction  of  scattering 
from  gaseous  bodies  with  high  aspect  ratios.  In  this  study,  the  FMM  is  directly  compared  with 
the  T-matrix  method  and  the  exact  prolate  spheroidal  solution  over  a  range  of  aspect  ratios 
up  to  20:1  for  constant  volume  prolate  spheroids.  Predictions  from  the  T-matrix  method  and 
FMM  are  compared  in  Fig.  13  for  aspect  ratios  up  to  16:1;  however,  results  for  aspect  ratios 
up  to  20:1  were  presented  by  Feuillade  and  Werby  (1994)  using  the  T-matrix  method,  as  well 
as  Ye  and  Hoskinson  (1998)  using  the  exact  prolate  spheroidal  solution.  In  addition,  low  ka 
resonance  scattering  by  gaseous  bodies  of  higher  aspect  ratio  (up  to  100:1)  was  calculated  with  no 
difficulties  relative  to  the  numerical  issues  cited  earlier.  These  ratios  are  much  higher  than  what 
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has  previously  appeared  in  the  literature  (Strasberg,  1953;  Weston,  1967;  Feuillade  and  Werby, 
1994;  Ye  and  Hoskinson,  1998).  The  FMM  accurately  predicts  the  scattering  in  the  resonance 
region — the  increase  in  the  resonance  frequency,  the  decrease  in  the  resonance  amplitude  and 
the  broadening  of  the  resonance  peaks  for  prolate  spheroids  of  constant  volume  and  increasing 
aspect  ratio  within  the  verifiable  range  of  aspect  ratios  presented  in  the  literature  (up  to  20:1 
aspect  ratios). 

V.  SUMMARY  AND  CONCLUSIONS 

An  extension  of  a  two-dimensional  conformal  mapping  approach  to  scattering  by  irregular, 
finite-length  bodies  of  revolution  has  been  presented  for  three  boundary  conditions.  The  model 
conformally  maps  the  coordinate  variables  of  the  original  coordinate  system  to  a  new  orthogonal 
coordinate  system  in  which  the  locus  of  all  points  where  the  new  radial  coordinate  is  a  con¬ 
stant  exactly  coincides  with  the  scatterer  surface.  The  solutions  to  the  transformed  Helmholtz 
equation  are  a  general  solution  for  the  total  pressure  in  the  case  of  far-field  scattering  from  a 
finite  body  of  revolution.  This  model  has  been  shown  to  be  very  accurate  in  the  prediction  of 
scattering  from  smooth,  symmetric  bodies  for  a  wide  range  of  frequencies  (Rayleigh  through  geo-  o* 
metric  scattering  region),  scattering  angles  (monostatic  and  bistatic),  aspect  ratios  and  boundary 
conditions,  and  for  all  angles.  Reasonable  agreement  has  also  been  demonstrated  for  irregular, 
realistic  shapes  when  compared  to  the  Kirchhoff  approximation.  This  work  represents  a  sig¬ 
nificant  advancement  by  providing  a  numerically  efficient  formulation  that  is  applicable  over  a 
wide  range  of  frequencies,  scattering  angles  and  geometries  for  soft,  rigid  and  fluid  boundary 
conditions. 

There  is  great  potential  for  further  development  of  this  model,  including  application  to  acous- 
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tic  scattering  by  elastic  solids  and  shells.  The  implementation  of  better  numerical  integration 
techniques,  the  use  of  extended  precision  format  of  floating  point  arithmetic  or  various  smooth¬ 
ing  techniques  (Yamashita,  1990)  could  delay  the  onset  of  ill-conditioned  matrices,  accelerating 
convergence.  New  scaling  techniques  based  upon  the  physical  scattering  mechanisms  could  also 
increase  accuracy  (Schmidt,  1993). 
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Figure  Captions 


Fig.  1:  Scattering  geometry  for  an  irregular,  axisymmetric  finite-length  body.  The  body 
is  symmetric  about  the  z-axis.  The  azimuthal  angular  coordinates,  <p  and  v,  range  from  0 
to  2  n  in  the  Ay-plane,  and  the  polar  angular  coordinates,  6  and  w,  range  from  0  to  it, 
measured  from  the  z-axis.  The  radial  coordinate  in  the  («,  w,  v)  coordinate  system  equals 
zero  on  the  surface.  Broadside  incidence  corresponds  to  6  -90  degrees.  End-on 
incidence  corresponds  to  0  -0  and  180  degrees.  In  the  new  coordinate  system,  g(u,w)  is 
the  length  along  the  z-axis,  mdj(u,w)  is  the  projection  in  the  Ay-plane. 

Fig.  2:  Performance  envelope  for  broadside  backscatter  for  a  smooth  prolate  spheroid 
with  soft  boundary  conditions  as  a  function  of  ka  and  aspect  ratio.  A  converged  solution 
is  defined  in  this  paper  as  one  in  which  the  computation  of  additional  modes  changed  the 
scattering  amplitude  by  less  than  0.1%  for  a  given  value  of  ka.  As  the  aspect  ratio  or 
value  of  ka  is  increased,  converged  solutions  are  more  difficult  to  obtain.  Aspect  ratio 
(AR)  is  the  ratio  of  length  to  width  of  the  prolate  spheroid  (AR=1  for  a  sphere).  The 
semi-minor  axis  of  the  prolate  spheroid  is  designated  as  a. 

Fig.  3:  Performance  envelope  for  broadside  backscatter  by  a  smooth  prolate  spheroid 
with  rigid  boundary  conditions  as  a  function  of  ka  and  aspect  ratio.  Converged  solutions 
are  defined  in  the  caption  of  Fig.  2. 


Fig.  4:  FMM  and  exact  solutions:  Reduced  target  strength  as  a  function  of  ka  for  soft, 
rigid  and  fluid  spheres.  Mass  and  sound  speed  contrasts  for  the  (weakly  scattering)  fluid 
case  arc  g=1.043  and  h=1.052,  respectively.  Exact  solution  calculation  based  on 
formulation  from  Anderson  (1950). 

Fig.  5:  FMM  and  DFC:  Reduced  target  strength  of  rigid  prolate  spheroids  at  broadside 
incidence  as  a  function  of  ka  for  aspect  ratios  of  2:1, 5:1  and  10:1.  DFC  calculations  are 
based  on  the  formulation  from  Stanton  (1989a).  The  agreement  between  the  FMM  and 
DFC  improves  as  the  aspect  ratio  increases,  which  is  consistent  with  the  fact  that  the  DFC 
is  valid  in  the  limit  of  high  aspect  ratio  geometries. 

Fig.  6:  FMM  and  T-matrix:  Magnitude  of  backscattered  and  bistatic  form  function  of  a 
2: 1  rigid  prolate  spheroid  as  a  function  of  ka  for  the  case  of  broadside  backscatter, 
backscatter  at  oblique  incidence  and  bistatic  scattering.  T-matrix  results  are  from  Figs.  3, 
4  and  5  of  Varadan  et  al.  (1982).  The  form  function  is  given  in  Eq.  (46). 

Fig.  7:  FMM  and  T-matrix:  Magnitude  of  backscattered  form  function  for  a  10:1  rigid 
prolate  spheroid  as  a  function  of  ka  for  the  case  of  broadside  incidence  and  end-on 
incidence.  T-matrix  results  are  from  Fig.  37(a)  of  Hackman  (1993).  The  form  function  is 
given  in  Eq.  (46). 


Fig.  8:  FMM  and  BEM:  Reduced  target  strength  as  a  function  of  ka  for  a  2:1  rigid 
prolate  spheroid  at  end-on  incidence.  The  agreement  continues  well  past  ka-  10.  BEM 
results  are  from  Francis  (2001). 

Fig.  9:  FMM  and  exact  prolate  spheroidal  solution:  Reduced  target  strength  as  a 
function  of  ka  for  a  10: 1  rigid  prolate  spheroid  at  incidence  angles  of  30, 70  and  90 
degrees  for  backscatter  and  forward  scatter.  Exact  prolate  spheroidal  solutions  are  from 
Fig.  7  of  Ye  et  al.  (1997). 

Fig.  10:  FMM  and  Kirchhoff  approximation:  Reduced  target  strength  as  a  function  of 
angle  for  a  5:1  soft  prolate  spheroid  at  ka  =1,  ka  =2  and  ka  =3.  Kirchhoff  approximation 
results  are  based  on  the  formulation  from  Bom  and  Wolf  (1999).  The  largest  differences 
occur  well  off  broadside,  where  the  Kirchhoff  approximation  characteristically 
underpredicts  the  scattering. 

Fig.  11:  FMM  and  Kirchhoff  approximation:  Target  strength  as  a  function  of  angle  for 
irregular  Shape  #1  at  ka  =1,  ka  =2,  ka  =3  and  ka  =4.  Kirchhoff  approximation  results  are 
based  on  the  formulation  from  Clay  and  Home  (1994).  Shape  #1  is  shown  in  the  bottom 
frame.  The  vertical  axis  is  exaggerated  to  better  illustrate  the  irregularity  of  the  surface. 
Broadside  corresponds  to  zero  degrees. 

Fig.  12:  FMM  and  Kirchhoff  approximation:  Target  strength  as  a  function  of  angle  for 
irregular  Shape  #2  at  ka  =1,  ka  =2,  ka=3,ka=4  and  ka  =5.  Kirchhoff  approximation 


results  are  based  on  formulation  from  Clay  and  Home  (1994).  Shape  #2  is  shown  in  the 
bottom  frame.  The  vertical  axis  is  exaggerated  to  better  illustrate  the  irregularity  of  the 
surface.  Broadside  corresponds  to  zero  degrees. 

Fig.  13:  FMM  and  T-matrix:  Magnitude  of  normalized  backscattering  amplitude  as  a 
function  of  kac  for  a  gaseous  (£=0.00126,  h=0.22)  prolate  spheroid  at  broadside 
incidence.  Results  from  the  FMM  (top  panel)  are  compared  to  results  from  the  T-matrix 
method  (bottom  panel)  from  Fig.  4(a)  of  Feuillade  and  Werbe  (1994).  The  single  dots  in 
the  top  panel  are  peak  amplitudes  for  aspect  ratios  listed  in  parentheses.  The  scattering 
amplitude  on  the  vertical  axis  is  normalized  to  the  peak  amplitude  of  the  sphere.  As  the 
aspect  ratio  of  the  spheroid  is  increased  beyond  1:1,  the  dimensions  of  the  spheroid  are 
changed  to  maintain  a  volume  equal  to  that  of  the  original  sphere.  Thus,  the  a*  in  ko*  is 
the  equivalent  spherical  radius  of  the  object,  corresponding  to  the  radius  of  the  sphere 
(the  1:1  case  shown)  that  has  the  same  volume  as  the  prolate  spheroid.  Results  for  aspect 
ratios  up  to  20:1  were  tabulated  by  Feuillade  and  Werby  (1994)  using  the  T-matrix 
method,  as  well  as  Ye  and  Hoskinson  (1998)  using  the  exact  prolate  spheroidal  solution. 
There  is  excellent  agreement  among  predictions  from  all  three  formulations  for  aspect 
ratios  up  to  20:1. 
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